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Abstract: We study QCD with two flavors of non-perturbatively improved Wilson 
fermions at finite temperature on the 16 3 8 lattice. We determine the transition 
temperature at lattice spacings as small as a ~ 0.12 fm, and study string breaking 
below the finite temperature transition. We find that the static potential can be 
fitted by a two-state ansatz, including a string state and a two-meson state. We 
investigate the role of Abelian monopoles at finite temperature. 
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1. Introduction 



Recently, efforts have been made to determine the critical temperature T c of the 
finite temperature transition in full QCD with Nf = 2 flavors of dynamical quarks. 
The Bielefeld group employed improved staggered fermions and an improved gauge 
field action [1]. The CP-PACS collaboration used improved Wilson fermions with 
mean field improved clover coefficient and an improved gauge field action [2]. Both 
groups were able to estimate T c in the chiral limit, and their values are in good 
agreement with each other. Still, there are many sources of systematic uncertainties. 
The main one is that the lattices used so far are rather coarse. In this paper we 
perform simulations on N% N t = 16 3 8 lattices at lattice spacings a much smaller 
than in previous works [1, 2]. To further reduce finite cut-off effects, we use non- 
perturbatively 0(a) improved Wilson fermions. A small lattice spacing is particularly 
helpful in determining the parameters of the static potential. 

In the presence of dynamical quarks the flux tube formed between static quark- 
antiquark pairs is expected to break at large distances. At zero temperature T 
the search for string breaking, i.e. flattening of the static potential, has not been 
successful, if the static quarks are created by the Wilson loop. (See e.g. [3]). At 
finite temperature string breaking has been observed at T < T c , if Polyakov loops 
are used instead to create the quarks. It is important to know the static potential 
at finite temperature for phenomenology [4]. In particular, it is needed to compute 
the dissociation temperatures for heavy quarkonia. We suggest a new ansatz and 
confront that ansatz with our numerical data. 

The dynamics of the QCD vacuum, and color confinement in particular, becomes 
more transparent in the maximally Abelian gauge (MAG) [5, 6]. In this gauge the 
relevant degrees of freedom are color electric charges, color magnetic monopoles, 
'photons' and 'gluons' [7]. There is evidence that the monopoles condense in the low 
temperature phase of the theory [6, 8], causing a dual Meissner effect, which constricts 
the color electric field into flux tubes, in accord with the dual superconductor picture 
of confinement [9]. Abelian dominance [10] and the dynamics of monopoles have been 
studied in detail at zero temperature in quenched [11] and unquenched [12, 13] lattice 
simulations. It turns out that in MAG the string tension is accounted for almost 
entirely by the monopole part of the Abelian projected gauge field [14, 15, 16, 12]. 
Furthermore, in studies of SU(2) gauge theory at nonzero temperature [17] it has been 
found that at the phase transition the Abelian Polyakov loop shows qualitatively the 
same behavior as the non-Abstain one. In this paper we extend the investigation of 
Abelian dominance to full QCD at nonzero temperature. 

The paper is organized as follows. In Section 2 we present the details of our 
simulation. Furthermore, we describe the gauge fixing algorithm and the Abelian 
projection. Section 3 is devoted to the determination of the transition temperature, 
and in Section 4 our results for the heavy quark potential are presented. In Section 
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5 we study the monopole density in the vacuum, as well as the action density in 
the vicinity of the flux tube. We demonstrate that the flux tube disappears at large 
quark- ant iquark separations. Finally, in Section 6 we conclude. Preliminary results 
of this work have been reported in [18]. 

2. Simulation details 

We consider Nf — 2 flavors of degenerate quarks. We use the Wilson gauge field 
action and non-perturbatively 0(a) improved Wilson fermions [19] 

S F = SP - g c sw a 5 ^2i!{s)a^F f , u (s)ilj{s), (2.1) 

s 

where is the original Wilson action, g is the gauge coupling and F fll/ (x) is the 
field strength tensor. The clover coefficient c sw is determined non-perturbatively. 
This action has been used in simulations of full QCD at zero temperature by the 
QCDSF and UKQCD collaborations [20, 21], whose results we use to fix the physical 
scale and the m^/mp ratio. At finite temperature the same action was used before 
in simulations on N t = 4 and 6 lattices at rather large quark masses and lattice 
spacings [22]. 

Non-perturbatively improved Nf = 2 Wilson fermions should not be employed 
below (3 = Q/g 2 = 5.2. In fact, c sw is known only for f3 > 5.2 [23]. The simulations 
are done on 16 3 8 lattices at two values of the coupling constant, (3 = 5.2 and 5.25, 
and nine different k values each. The parameters are listed in Table 1. They are also 
shown in Fig. 1, together with lines of constant r Q /a and constant m n /m p obtained 
at T = [20]. Note that the lines of constant T run parallel to the lines of constant 
r /a. To check the finite size effects we have also done simulations on the 24 3 8 lattice 
at = 5.2, k = 0.1343. 

The dynamical gauge field configurations are generated on the Hitachi SR8000 at 
KEK (Tsukuba) and on the MVS 1000M at Joint Supercomputer Center (Moscow), 
using a Hybrid Monte Carlo , while the analysis is done on the NEC SX5 at RCNP 
(Osaka) and on the PC-cluster at ITP (Kanazawa). Our present statistics is shown 
in Table I. The length of the trajectory was chosen to be r = 0.25. We use a blocked 
jackknife method to compute the statistical errors of the observables and a bootstrap 
method to compute the errors of the fit parameters. We compute the Polyakov loop 

1 Nt 

^ = 3^II^ 4 )' ( 2 - 2 ) 

S4 = l 

U (s, ji) being the link variable, on every trajectory. From that we derive the suscep- 
tibility 

x = N*j2(( L2 ^)-(mn (2.3) 
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Figure 1: Lines of constant ro/a (solid lines) and constant m^/mp (dotted lines) at T = 0. 
Crosses indicate points where simulations are done. 
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46 


0.1350 


1780 


30 



Table 1: Parameters and statistics of the simulation, together with the integrated auto- 
correlation time. The length of the trajectory is r = 0.25. 

and the integrated autocorrelation time Ti nt . The autocorrelation time is given in Ta- 
ble 1 in units of trajectories. Furthermore, we compute the Polyakov loop correlator 
{L(s)L^(s ')), from which we obtain the static potential. To reduce the error on the 
static potential, we employ a hypercubic blocking of the gauge field as described in 
[24]. We choose every 5 th to 20 th trajectory, depending on the value of k, to compute 
the blocked Polyakov loop correlator. 
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We fix the MAG [6] by maximizing the gauge fixing functional F(U), 

F ( f/ ) = I^7E(l^( S ^)| 2 +l^2( S ,/i)r+|f/33( S ,/i)| 2 ) (2.4) 

with respect to local gauge transformations g of the lattice gauge field, 

U(s, p) - U%s, p) = g( S yU(s, p)g{s + p) . (2.5) 

To do so, we use the simulated annealing (SA) algorithm. The advantage of 
this algorithm over the usual iterative procedure has been demonstrated in [16] for 
the MAG and in [25] for the maximal center gauge in the SU(2) gauge theory. In 
practice one does not find the global maximum of the gauge fixing functional in a 
finite-time computation. For this reason it has been proposed [16] to apply the SA 
algorithm to a number of randomly generated gauge copies and pick that one with 
the largest value of F. By increasing the number of gauge copies, one eventually 
reaches the situation, where the statistical noise is larger than the deviation from 
the global maximum. We use one gauge copy at (5 = 5.2 and three gauge copies at 
f3 = 5.25. We have checked that by increasing the number of gauge copies our results 
for the gauge dependent quantities are left unchanged within the error bars. 

To obtain Abelian observables, one needs to project the SU(3) link matrices onto 
the maximal Abelian subgroup U(l) x U(l) first. The original construction [26] is 
equivalent to finding the Abelian gauge field u(s,p) £ (7(1) x U(l) which maximizes 
|Tr (U(s, fi)vt(s, p)) | 2 . The Abelian counterpart of an observable is then obtained 
by substituting 

u(s, p) = diag(e Wl{s ^\e W2{s ^\e W3{s '^) , (2.6) 

for U(s,p); X!i=i ^i( s > A 4 ) = 0, so that det(u(s, p)) = 1. From eq. (2.6) we define 
plaquette angles 

0i(s, pv) = dpOifa v) - dJi(s, p) , (2.7) 

where <9 M is the lattice forward derivative. The plaquette angles can be decomposed 
into regular and singular components, 

9i(s, \w) = 9i(s, \w) + 2irnii(s, \xv) , (2.8) 

where 9i(s,fiu) e (— 7r;7r], and mi(s,fiu) e N counts the number of Dirac strings 
piercing the given plaquette. Note that Yli 9i( s , = 2tt/, I — 0, ±1. If I — +1 (—1) 
we substitute the largest (smallest) 9i(s, pa/) (of the three components) by 9i(s, pa/) — 
2n (+27r), and similarly for rrii{s,pv), so that J2i9i(s, pu) = J2i m i( s i A^) = 0- 

The monopole currents, being located on the links of the dual lattice, are defined 
by [27] 

ki(*s, p) = —e^ pa d v 9i(s + p, pa) = -\t^ pa d u mi(s + p, pa) . (2.9) 
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They satisfy the constraint 



^r^(* s ,/i) = o, (2.10) 

i 

for any *s,/i. The Abelian gauge fields 9i(s,fi) can in turn be decomposed into 
monopole (singular) and photon (regular) parts: 

ofaii) = ri + 0?(s,ri- (2-ii) 

The monopole part is defined by [28]: 

er n (s, fi) = -27T D ( s - s ') K Ms', "I*) , (2.12) 

s' 

where d v ' is the backward lattice derivative, and D(s) denotes the lattice Coulomb 
propagator. 

The Abelian Polyakov loop is defined by 

1 3 N t 

LaU?) = gE^ 5 )' ^f bcl ( S 1 =exp{i^^( S ,4)}. (2.13) 

i=l S4=l 

Similarly, we define the monopole and photon Polyakov loops [29]: 

, 3 Nt 

= 3 J2 L ™ a W » L ™ n & = ^ E d r\sA)} , (2.14) 
1=1 S4=l 

V(^) = l E ^(^) , ^f(^) = exp{i E ^(^ 4)} . (2.15) 

i=l S4=l 



3. Transition temperature 

The order parameter of the finite temperature phase transition in quenched QCD is 
the Polyakov loop, and the corresponding symmetry is global Z(3). In the presence 
of dynamical 'crural' fermions the chiral condensate (ipip) is an order parameter (of 
the chiral symmetry breaking transition). It is expected that there is no phase 
transition at intermediate quark masses, only a crossover. Numerical results show 
that both order parameters can be used to locate the transition point at intermediate 
quark masses [1]. We use the Polyakov loop, because the calculation of the chiral 
condensate for Wilson fermions requires renormalization and is rather involved. 

It is instructive to plot the Polyakov loop in the complex plane as a function of 
temperature, which has been done in Fig. 2. We find that the distribution is rather 
asymmetric, even at the lowest temperature, favoring a positive value of Re L. This 
is indeed what one expects [30]. The introduction of dynamical quarks adds a term 
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Figure 2: Scatter plots of the Polyakov loop in the complex plane for various temperatures. 
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Table 2: The expectation values of the non-Abelian, Abelian, monopole and photon 
Polyakov loops at (3 = 5.2. 

proportional to ReL to the effective action, which results in a nonzero value of (L). 
The numbers are given in Tables 2 and 3. In the Tables we also give values for the 
Abelian and monopole Polyakov loops separately. 

As we can see from Fig. 1, increasing ft at a fixed value of (3 increases the 
temperature T oc r /a. In Figs. 3, 4, 5 we plot the expectation values of the various 
Polyakov loops of Tables 2 and 3 as a function of k. While (L), (LAbei) and (L mon ) 
increase with increasing k, (L p h) stays approximately constant over the full range of k. 
Furthermore, similar to the quenched theory, L mon and L ph are virtually independent, 
which follows from (LAbei) ~ (L moTL ) (£ p h) • We have no explanation for the dip seen 
in (L) at (5 — 5.25 around k = 0.1341. We see some signal of metastability, but 
we will need higher statistics to clarify this point. A similar dip is seen on Fig. 1 
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Figure 3: The expectation value of the non-Abelian Polyakov loop as a function of k for 
(3 = 5.2 and 5.25. 
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Table 3: The expectation values of the non-Abelian and monopole Polyakov loops at 
(3 = 5.25. 



of ref. [22], where the same lattice action was used to study the phase transition on 
small lattices. 

The task is now to determine the transition temperature T c . We call the k 
value, at which the transition takes place, Kt- We identify K t as the point, where the 
Polyakov loop susceptibility (2.3) assumes its maximum. The Abelian, monopole 
and photon Polyakov susceptibilities XAbei, Xmon and Xph, respectively, are defined 
similarly. The susceptibilities are given in Tables 4 and 5, and they are plotted in 
Figs. 6 and 7. From the non-Abelian susceptibility \ we find n t = 0.1344(1) at 
(3 = 5.2 and n t = 0.1341(1) at (3 — 5.25 where the central values and the errors 
are determined by the maxima of susceptibilities and by the distances to neighbor 



-7- 



<L> 
0.6 

0.4 

0.2h 





<$> Monopole, [3=5.25 
□ Monopole, p=5.2 
- < Photon, [3=5.2 $ 



-| 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 r 

□ □ 



* ID m 



_i i i i I i i i i I i i i i I i i_ 



0.133 0.134 0.135 0.136 K 



Figure 4: The expectation value of monopole and photon Polyakov loops as functions of 
k for p = 5.2 and 5.25. 
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Figure 5: The expectation value of the Abelian Polyakov loop as a function of k for 
(3 = 5.2. The product of monopole and photon Polyakov loops is also shown. 



data points, respectively This translates into T c r$ = 0.53(1) at /3 — 5.2 and 
T c tq = 0.56(1) at (3 — 5.25, where the numbers of r /a at K t have been obtained by 
interpolation of the T = results [20]. Taking r = 0.5 fm to fix the scale, we obtain 
in physical units 

T c = 210(3) MeV, (3 = 5.2 (3.1) 
T c = 219(3) MeV, (3 = 5.25. 
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Figure 6: The non-Abelian Polyakov loop susceptibility as a function of k for (3 = 5.2 
and 5.25, respectively, together with the fit. 
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Table 4: The same as in Table 3 but for the susceptibility of the Polyakov loop at (3 = 5.2. 



By interpolating m^/rrip, given in [20], to n t we obtain at the transition point 
m^/rrip = 0.77 at (3 = 5.2 and 0.81 at (3 — 5.25. Similarly, we can compute the 
temperature T at our various k values. The result is given in the last column of 
Table 2 and Table 3 in the form of T/T c . 

The susceptibilities XAbei an d Xmon have maxima at the same k value as the 
non-Abelian susceptibility. We find that XAbei ~ (£ph) 2 Xmon- This follows from our 
earlier observation, namely that L mon and L pn are independent, and the smallness 
of Xph- The non-Abelian susceptibility is 10 to 50 times smaller than its Abelian 
counterpart. The photon susceptibility does not show any change at the critical 
temperature, as expected. We conclude, that the monopole degrees of freedom are 
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Figure 7: The same as in Fig. 6 but for the Abelian, the monopole and the photon 
Polyakov loop susceptibilities. The Abelian (photon) Polyakov loop susceptibility has been 
enhanced by a factor of 10 (100). 
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Table 5: The same as in Table 4 but for the susceptibility of the Polyakov loop at [3 = 5.25. 

most sensitive to the transition, as was the case in the quenched theory. 

In Fig. 8 we compare our results for T c with those of Refs. [1] and [22], where 
we have assumed \fd = 425 MeV. Our results are in quantitative agreement with 
the results of the Bielefeld group. This is reassuring, as [1] and [22] work at larger 
lattice spacing. 

We fit the susceptibility in the transition region by [31] 



(3.2) 



C 2 +(K-Kt) 2 

where C\ 2, Kt and v are taken as fit parameters. The fit values, and the corresponding 
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Figure 8: The transition temperature as a function of m^/mp from this work (triangles), 
and from [1] (squares) and [22] (circle), respectively. 
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Table 6: The parameters of the fit (3.2). The critical temperature T c , determined from 
fits of various susceptibilities, is given in MeV. 
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values for T c , are presented in Table 6. 

The fits are shown by the dashed lines in Fig. 6. All fits gave x 2 /dof = 0(1). 
We observe that the non-Abelian, Abelian and monopole susceptibilities (available 
at (3 = 5.2) give the same n t within error bars. In principle, this fit should give the 
same value of Kt as our previous estimate by the maximum of susceptibility but with 
better precision. This is indeed the case for (3 = 5.25. For (3 = 5.2 the result is 
shifted to higher k. This shift is due to the data at k — 0.1348, which is rather far 
from the transition. To be on a safe side we take the values quoted in eq.(3.1) as the 
critical temperature and respective errros for both values of (3. 

Because N s /N t = 2 in our simulations, the question of finite volume effects 
is essential. To check for finite volume effects we have performed simulations at 
(3 = 5.2, k = 0.1343 on 24 3 x 10 lattices. The value of k was chosen close to the 
transition point, where finite volume corrections are expected to be largest. We 
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found (L) = 0.0098(10) and x = 0.099(13), as compared to (L) = 0.0092(13) and 
X = 0.115(17), respectively, on the N s = 16 lattice. In both cases the numbers agree 
within the error bars, so that we do not reckon with large effects. 



4. Heavy quark potential 

4.1 Ansatz 

One of the characteristic features of full QCD is breaking of the string spanned 
between static quark and anti-quark pairs. String breaking will manifest itself in a 
type of screening behavior of the heavy quark potential. At zero temperature no 
clear evidence for string breaking has been found (in QCD) yet. The reason is, so it 
is believed, that the Wilson loop has very small overlap with the broken string state. 
The expectation value of the Wilson loop for large distances r can be written as 

(W(r,t)) = C v (r)e' {Vo+v ^^ r))t + C E (r)e- 2E ^ + ... } (4.1) 

where V^ tr ing( r ) is the usual confining potential, V^ tr ing(^) — —a/r + ar, E s \ is the 
static-light meson energy, and Vo is the self-energy. The latter is divergent in the 
continuum limit a — > 0. The overlap with the string state, Cy(r), is of the order 
of one, while the overlap with the broken string, C^(r), appears to be small. An 
estimate [32] is: C E (r) ~ e" 2msl ' r , where m s i = E si — V /2 is the so called binding 
energy of the static-light meson or, in other words, the constituent quark mass [33]. 
(See also the discussion in Ref. [34].) A similar estimate was given in Ref. [35], based 
on the hypothesis of Abelian dominance. 

The conventional definition of the string breaking distance is the distance r s t, at 
which the energy of two static-light mesons is equal to the energy of the string, i.e. 

71 

2m s i — (7 • r s b — — . (4.2) 

The Wuppertal group found r sb = 2.3 r at m^/rrip = 0.7 [3], while CP-PACS found 
r s b = 2.2 r at m^/rrip = 0.6 [32]. In full QCD it was found y/ar = 1.14 [3] and 
= 1.16 [21], respectively, from which we derive 2m s i m 2.9/r ~ 1.1 GeV, 
assuming tq = 0.5 fm. This agrees with the estimate of [36]. Using these values for 
a/cxto and m s i and assuming the mentioned above form of Ce{t) we can estimate the 
values of r and t at which the two terms in eq.(4.1) become equal indicating that the 
string breaking effects become visible: 

r = t = ^. (4.3) 
o 

We also estimate the numerical value of the Wilson loop of the corresponding size: 

(W(r,t)) < 10- 11 -e- Vot . (4.4) 
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It is a challenging task to record a Wilson loop of this order of magnitude. Recently, 
a first successful attempt to do so was reported in [37] where the authors studied the 
adjoint static potential in three-dimensional SU(2) gauge theory. 

At finite temperature T < T c string breaking has been studied in [38] . The heavy 
quark potential V(r,T) is obtained from the Polyakov loop correlator: 

±V(r,T) = -\n(L(s)L\s')), (4.5) 
up to an entropy contribution, where r = \s — s'\. At large separations 

(mtf(s')) j \(L)\\ (4.6) 

where \(L)\ 2 ^ 0, as global Z z is broken by fermions. It should be noted that the 
potential in (4.5) is a color average, which is related to the proper singlet and the 
octet potentials by [39, 40] 

e -V(r,T)/T = \ e -V sing (r,T)/T + § e -V oct (r,T)/T _ ^ ^ 

It would be desirable to compute the singlet potential introduced in (4.7). Work on 
calculating the singlet and octet potential separately is in progress. In the recent 
publication [41] this calculation was already performed for Nf = 2 lattice QCD on 
16 3 4 lattice. 

The spectral representation of the Polyakov loop correlator is given by [42] 

oo 

(L(s)L\s')) = J2^ne- E - {r)/T , (4.8) 

?1=0 

where w n are integers. At zero temperature we have V(r, 0) = E (r), up to a 
constant, where E is the ground state energy. At finite temperature V(r,T) gets 
contributions from all (excited) states. As was discussed already, at T = the 
potential can be described by the string model potential up to the string breaking 
distance r s b. Beyond this distance the state of two static-light mesons becomes the 
ground state of the system. Thus, there are two competing states in the spectrum, 
and it depends on the distance r, which one will be the ground state. We may expect 
that the situation at small temperature is similar to the case of T = 0. 

We shall now assume that at temperatures T < T c the Polyakov loop correlator 
can be described in terms of these two states. We then have 

(L(s)L\s')) = e -(Vo(T)+V string (r,T))/T + e -2E(T)/T ^ (4Q) 

where the finite temperature string potential V r str i ng (r, T) is given by [43]: 

14tring(r,T) = --(a- iarctan(2rT)) 

/ 2T 2 1 \ T / \ 

+ (a(T) + — arctan^Jr + - \n(l + Ar 2 T 2 j . (4.10) 
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We consider the temperature dependent string tension cx(T) as a free parameter. 
While in [43] a was fixed at tt/12, a fit of the short distance part of the potential 
at T = gave [44] a = 0.32 ~ 0.34. In the following we shall consider both cases, 
a = 7r/12 and a = 0.33. The energy E(T) can be written as 

E(T) = ±V + m(T), (4.11) 

where m(T) is the constituent quark mass [33]. 

A long time ago the following ansatz for the Polyakov loop correlator has been 
proposed [4]: 

(L(s)L^s')) = e-( y o(T)+vkM S (r,T))/T ? (412) 

where 

V KMS (r,T) = -(1 - - -e-» r . (4.13) 
\i r 

As we will see, this potential cannot capture the physics below the finite temperature 
transition. We also do not consider this potential a valid ansatz, because string 
breaking is a level crossing phenomenon. 

Besides the non-Abelian potential, we will study the Abelian one. In particu- 
lar we shall be interested in its monopole and photon parts. From studies at zero 
temperature [16, 12] it is known that the monopole part of the potential decreases 
linearly down to very small distances, showing no Coulomb term, which sometimes 
makes it easier to extract a string tension. It appears that the monopole part of the 
potential has not only no Coulomb term, but also shows no broadening of the flux 
tube as the length of the flux tube is increased. (Both phenomena are connected of 
course.) As we show below, our monopole potential is also linear at distances up to 
the distance of order of 0.5 fm where flattening starts. Thus we may write 

< W^Ln(s')) = e-<™^ M/T + e-* E -^ T , (4.14) 

where 

KSO, T) = a mon ■ r , E mon (T) = V mon (T) + m mon (T) . (4.15) 
4.2 Hypercubic blocking 

As was mentioned already, we apply hypercubic blocking (HCB) [24] to reduce the 
statistical errors. That means every SU(3) link matrix U(s,fi) is replaced by a new 
link matrix Uhcb(s, aO, which is the weighted sum of products of link matrices along 
paths from s to s + ft within adjacent cubes projected onto the nearest SU(3) group 
element. We used the same parameters as in [24]. 

In Fig. 9 we compare the static potential from blocked and unblocked configu- 
rations. We see that the statistical errors are substantially reduced. Furthermore, 
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Figure 9: Effect of hypercubic blocking on the potential at f3 = 5.2, k = 0.1335. 

rotational invariance is improved, in agreement with earlier observations [24]. The 
blocking procedure decreases the self energy of the static sources, which causes a 
constant shift in the potential. In Fig. 9 we shift the potential by 1.33/r , so that 
it agrees with the unblocked potential at r = y/2a. We find good agreement at all 
distances, except perhaps at r = a. The shift agrees with the change in the asymp- 
totic value of the potential, — 2Tlog(L), which was found to be 1.24(25)/r . The 
discrepancy at r = a can be accounted for by perturbative corrections [46]. All our 
fits are made for r/ro > 1 thus this point is always discarded. 

4.3 Non-Abelian potential 

We first fit the static potential by the two-state ansatz (4.9). This is done for 
two different choices of a, a = it/ '12 and a = 0.33. Examples of the fit for Tj 
T c = 0.87 and 0.98 and the second choice a = 0.33 are shown in Fig. 10. The 
curves for a = n/12 and a = 0.33 are practically indistinguishable from each other, 
visually and in terms of x 2 /dof. We also show the asymptotic value of the potential, 
— 2Tln(L). The potential converges to this value at large distances. The two-state 
ansatz describes the data very well. The fit parameters are given in Tables 7 and 8, 
where cr(0) = (1.14/r ) 2 [3] was used. 

The string tension cx(T) and the constituent quark mass m(T) are plotted in 
Figs. 11 and 12, respectively. In Fig. 11 the quenched value of a(T)/a(0) [47] is shown 
for comparison. Both the string tension and the constituent quark mass decrease 
with increasing temperature, as we expect. The results differ by approximately one 
a between a = 0.33 and a = n/12. For lower temperatures, at T/T c = 0.80 and 
T/T c = 0.87, we find rather good agreement between the results of the quenched 
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Figure 10: The heavy quark potential at [3 = 5.2 for T/T c = 0.87 and 0.98, together with 
the fit using a = 0.33. The horizontal lines show the asymptotic value of the potential, 
where the shaded area indicates the error. 



a = 0.33 


T/T c 


V r 


a(T)/a(0) 


m(T)r 


0.80 


2.33(2) 


0.86(3) 


1.02(10) 


0.87 


2.59(3) 


0.76(4) 


0.89(10) 


0.94 


2.80(3) 


0.79(5) 


0.70(8) 


0.98 


2.85(8) 


0.90(12) 


0.25(6) 




a = 7i 


/12 




0.80 


2.00(2) 


0.97(3) 


1.19(11) 


0.87 


2.23(3) 


0.89(4) 


1.00(10) 


0.94 


2.42(4) 


0.95(5) 


0.90(8) 


0.98 


2.46(10) 


1.05(14) 


0.44(7) 


Monopole part 


0.80 


0.47(1) 


0.90(1) 


1.24(4) 


0.87 


0.53(1) 


0.85(1) 


1.07(2) 


0.94 


0.61(1) 


0.84(1) 


1.03(3) 


0.98 


1.06(1) 


0.46(1) 


0.21(1) 



Table 7: Fit parameters of the two-state ansatz (4.9) at j3 
a(0) = (1.14/V ) 2 . 



5.2, where we have assumed 



theory and our results obtained with the choice a = 0.33. For higher temperatures 
agreement is much worse, especially for the lighter quark mass. A reason for this 
discrepancy might be that our Ansatz is not valid at temperatures close to the 
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a = 0.33 


T/T c 


V r 


a(T)/a(0) 


m(T)r 


0.86 


2.71(3) 


0.74(5) 


0.81(7) 


0.92 


2.87(7) 


0.70(10) 


0.37(7) 


0.96 


3.00(6) 


0.67(9) 


0.25(7) 


0.97 


3.12(7) 


0.57(10) 


0.15(7) 


0.99 


3.16(5) 


0.55(8) 


0.28(7) 


a = 7r/12 


0.86 


2.34(3) 


0.88(5) 


1.00(6) 


0.92 


2.49(7) 


0.86(10) 


0.57(7) 


0.96 


2.61(6) 


0.83(10) 


0.45(7) 


0.97 


2.71(6) 


0.73(9) 


0.35(6) 


0.99 


2.75(5) 


0.72(8) 


0.48(7) 


Monopole part 


0.86 


0.56(2) 


0.80(3) 


1.18(8) 


0.92 


0.57(1) 


0.80(2) 


0.54(2) 


0.96 


0.64(1) 


0.78(2) 


0.33(1) 


0.97 


0.66(1) 


0.79(1) 


0.27(1) 


0.99 


0.65(1) 


0.79(1) 


0.39(1) 



Table 8: Fit parameters of the two-state ansatz (4.9) at (3 = 5.25, where we have assumed 
a(0) = (1.14/ro) 2 . 

transition. 

Our values for the constituent quark mass m(T) are larger by about 100 MeV 
than those reported in [36]. However, one should note a difference between our 
definition of the self energy and the definition used in Ref. [36]. 

We are now able to compute the string breaking distance r sb from 

V stling (r sb ,T) = 2m(T). (4.16) 

In Figs. 13 and 14 we show the two energy levels together with the data. The string 
breaks where the two levels cross. The dependence of r sb on the temperature is 
shown in Fig. 15. We see that r sb decreases as the temperature is increased. The 
difference of r sb between the two choices of a lies within the error bars. 

Let us now consider the screening potential (4.13). Fitting this potential to 
our data gives a comparable value of x 2 /dof. However, the parameters of the fit 
turn out to be unphysical. For example, at T/T c = 0.80, 0.87 and 0.94 we obtain 
cr/a-(0) = 21(6), 13(2) and 5.5(6), respectively. Only close to the deconfmement 
transition we do find a reasonable value for the string tension: a/a(0) = 0.4(3) at 
T/T c = 0.98. 
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Figure 11: The string tension from a fit of the two-state ansatz (4.9) as a function of 
temperature. The quenched value of the string tension [47] is shown for comparison. The 
shaded region indicate the error bar. 
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Figure 12: The constituent quark mass from a fit of the two-state ansatz (4.9) as a 
function of temperature. The dash-dotted line indicate the zero-temperature quenched 
value of the mass. 



The screening potential (4.13) may be rewritten (up to a constant) in the fol- 
lowing form [45]: 

4 a, b(T)' 
3 r fi 



V Wong {r,T) 



-/tor 



(4.17) 



Taking (as in [45]) b(T) = b (1 - (T/T c ) 2 ), b = 0.35 GeV 2 , /i = 0.28 GeV, and 
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Figure 13: The string potential and the constituent mass as a function of distance at 
(3 = 5.2 for k = 0.1335 (T/T c = 0.87). The shaded regions indicate the errors. 
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Figure 14: The same as in Fig. 13 but for k = 0.1343 (T/T c = 0.98). 

a s ~ 0.32 (a s ~ 0.24) for the charmonium (bottonium) potential, and shifting the 
potential by a constant so that it agrees with the lattice potential at r = r , we find 
no agreement between this potential and the lattice data. Thus, the quarkonium 
spectra derived from this potential [45] need to be revised. 

4.4 Monopole part of the potential 

We carried out a similar analysis as before for the monopole part of the heavy quark 
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Figure 15: The string breaking distance function of temperature. The dash- 

dotted line indicate the corresponding zero-temperature quenched value. 

potential, which is obtained from the correlator (4.14). The fit parameters are given 
in Tables 7 and 8, and the potential is shown in Fig. 16. The errors are smaller than in 
the previous fits, as expected. The monopole part of the potential shows no Coulomb 
term, while at large distances it converges to its asymptotic value — 2Tln(L mon ). 
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Figure 16: The monopole part of the potential as a function of distance at (3 = 5.2 for 
T/T c = 0.87 and 0.98, together with a fit of the form (4.15) (dashed curve). The horizontal 
lines show the asymptotic value of the potential, where the shaded area indicates the error. 



The string tension cr mon (T) and the constituent quark mass m mon (T) are shown 
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Figure 17: The string tension of the monopole part of the potential as a function of 
temperature. In the left figure the quenched value of the string tension [47] is shown for 
comparison. The shaded region indicate the error bar. 
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Figure 18: The same as in Fig. 17 but for the constituent quark mass. The dash-dotted 
line indicate the zero-temperature quenched value of the mass. 



in Fig. 17 and 18, respectively. Because the Coulomb term is absent, we now can 
determine the string tension much more accurately. We find substantially larger 
values than in the quenched case. Furthermore, the string tension appears to decrease 
more slowly as the the system is heated. The constituent quark mass m mon (T) looks 
very much the same as in the non-Abelian case. The same holds for the string 
breaking distance, r^ on , which is shown in Fig. 19. 
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Figure 19: The string breaking distance r^ on obtained from the monopole potential 
as a function of temperature. The dash-dotted line indicate the corresponding zero- 
temperature quenched value. 



To shed further light on the string breaking mechanism, we have computed the 
action density, the color electric field and the monopole current in the vicinity of the 
(broken) string. 

The definitions of observables are the same as in Ref. [12]. We are interested in 
local Abelian operators of the form: 

0(s) = diag(Oi(s), 2 (s), 3 (s)) G U(l) x U(l) . (4.18) 

The correlator of the action density - which is C-parity even operator - with the 
product of the monopole Polyakov loops, Crao^s^C^^s"), can be written analo- 
gously to Ref. [48]: 

{ ° {S)) ^ = 3 (Tr[£ mon ( S 0/:L n ( S -)]) 3 (W) ' (419) 

where 

C mon (s) = dmg(Lr n (s),LT n (s), Lr n (s)) , (4.20) 

(cf., Eq. (2.14)). 

As for the C-parity odd operators O, such as the color electric field and the 
monopole current, we have 

_ WOW [A„o,(.?)£L„(?")])) ,, > 
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in analogy to the case of SU(2) and U(l) theories, Refs. [49]. 

The monopole part of the action density p^ mon , the monopole part of the color 
electric field Ef inon and the monopole current k Cin ° n , induced by the Polyakov loops, 
are then given by 

Pa—(s) = f $>iag(cos(£r n (s, I*, i/)), cos(^ on («, I*, ")), cos(£ 3 mon ( S , u)))) Cmon , 

ii>v 

(4.22) 

where the plaquette angles, 0™ on (s, /i, v), are constructed from the monopole link 
angles (2.12), 

Ef—{s) = i(diag(^r n (s,4,j)^ 2 mon (s,4,j)^3 mon ( S ,4,j)))£ mon , (4.23) 

and 

k c (*s, /i) = 2m (diag(A;i(*s, pi), k 2 (*s, pi), k 3 (*s, /i))) £mon , (4.24) 

respectively. 

In Fig. 20 we show the result for T/T c = 0.98 and three different separations, r = 
0.5, 0.8 and 1.3 fm. Our estimate of the string breaking distance at this temperature 
is 0.5 fm. The figure suggest that the flux tube has disappeared at the latest at 
r = 1.3 fm. 

5. Monopole density 

Another characteristic quantity of the confining vacuum is the monopole density, 
which we define as 

s i=l S,fM 

where the monopole current, fcj(*s, //), is given in (2.9). 

In Fig. 21 we compare the monopole density of this work with that of the 
quenched theory. The quenched result has been obtained on the same sized lat- 
tice at f3 — 5.8, 5.9, 6.0, 6.1 and 6.2. The density in full QCD is substantially higher 
than in the quenched theory, in agreement with our earlier result at T = [12]. We 
believe that the introduction of dynamical fermions causes an attraction between 
monopoles and antimonopoles, which naturally leads to an increase in the monopole 
density. A similar mechanism has been observed in the case of instantons and anti- 
instantons [50]. Both mechanisms are, of course, related, because (anti-) instantons 
are intimately connected with monopoles [51]. 
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Figure 20: The monopole part of the action density (a), the monopole part of the color 
electric field (b) and the solenoidal monopole current in the plane perpendicular to the flux 
tube (c) at T/T c = 0.98 and distances (from top to bottom) 0.5, 0.8 and 1.3 fm. 



Near the finite temperature transition we expect the monopoles to gradually 
become static as the temperature becomes high. This can be monitored by the 
asymmetry of the density of spatial and temporal monopole currents [26, 52]: 



V = 



Pt ~ Ps 



(5.2) 



Pt + Ps 

where pt(p s ) is the density of the temporal (spatial) monopole currents, 

3 3 3 



3 AW 



i=l s 



i=l s p,=l 



If all currents are time-like, then this quantity is unity, while in the case of an 
isotropic distribution it is zero. In Fig. 22 we plot the asymmetry r\ as a function 
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Figure 21: The monopole density at (3 = 5.2 as a function of temperature. 
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Figure 22: The asymmetry of the monopole density at (3 = 5.2 as a function of tempera- 
ture. 



of temperature. We compare the result with the predictions of the quenched theory. 
It is found that rj is zero in the confined phase and nonzero in the deconfined phase. 
In the deconfinement phase the value of rj is about 5 times smaller in full QCD 
compared to the quenched theory. A reason for this may be rooted in a different 
nature of the transition which is of the first order phase in the quenched case while 
in full QCD one observes a smooth crossover. 
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6. Conclusions 



We have studied QCD with two flavors of dynamical quarks at finite temperature 
on a 16 3 8 lattice. At the phase transition the lattice spacing is a ss 0.12 fm. We 
employed non-perturbatively improved Wilson fermions, so that we may expect finite 
cut-off effects to be small. 

To determine the parameters of the transition, notably the transition tempera- 
ture and the string tension, and to shed light on the dynamics of the transition, it 
helped to resort to Abelian variables in the maximally Abelian gauge. 

We observed string breaking in Polyakov loop correlators. This is a level crossing 
phenomenon. Accordingly, we fitted the correlator by a two-state ansatz, consisting 
of a string state and a two-meson state. We found good agreement of this ansatz with 
our numerical data for T < T c , while we could rule out previously proposed single- 
state correlation functions. The string breaking distance was found to be r sb ps 1 
fm at T/T c ps 0.8, our lowest temperature. String breaking is also clearly visible in 
the action density, the color electric field distribution and the solenoidal monopole 
current around the static sources. 

To make contact with the chiral limit, we need to increase N t , because the 
coupling cannot be taken smaller than (3 ps 5.2. Work on 24 3 10 lattices at m^/rrip Pi 
0.6 is in progress. Preliminary results have been presented in [53]. 

Acknowledgments 

We like to thank Alan Irving and Dirk Pleiter for assistance. The calculations have 
been performed on the Hitachi SR8000 at KEK Tsukuba and on the MVS 1000M 
at Moscow. We like to thank the staff of the Moscow Joint Supercomputer Center, 
especially A.V. Zabrodin, for their support. We furthermore thank Ph. De For- 
crand, V. Mitrjushkin and M. Muller-Preussker for useful discussions. This work is 
partially supported by grants INTAS-00-00111, RFBR 02-02-17308, RFBR 01-02- 
17456, RFBR-DFG 03-02-0416, RFBR 03-02-16941, DFG-RFBR 436RUS1 13/739/0 
and CRDF award RPI-2364-MO-02. M.Ch. is supported by JSPS Fellowship P01023. 
V.B. acknowledges support from JSPS RC30126103. T.S. is supported by JSPS 
Grant-in-Aid for Scientific Research on Priority Areas 13135210 and 15340073. 

References 

[1] F. Karsch, E. Laermann and A. Peikert, Nucl. Phys. B 605 (2001) 579. 

[2] A. Ali Khan et al. [CP-PACS Collaboration], Phys. Rev. D 63 (2001) 034502. 

[3] B. Bolder et al, Phys. Rev. D 63 (2001) 074504. 

[4] F. Karsch, M. T. Mehr and H. Satz, Z. Phys. C 37 (1988) 617. 



- 26 - 



[5] G. 't Hooft, Nucl. Phys. B190, 455 (1981). 

[6] A.S. Kronfeld, M.L. Laursen, G. Schierholz and U.-J. Wiese, Phys. Lett. B198, 516 
(1987). 

[7] A. S. Kronfeld, G. Schierholz and U.-J. Wiese, Nucl. Phys. B 293 (1987) 461. 

[8] H. Shiba and T. Suzuki, Phys. Lett. B 351 (1995) 519; M. N. Chernodub, 

M. I. Polikarpov, A. I. Veselov, Phys. Lett. B 399, 267 (1997); A. Di Giacomo and 
G. Paffuti, Phys. Rev. D 56, 6816 (1997). 

[9] G. 't Hooft, in High Energy Physics, ed. A. Zichichi, EPS International Conference, 
Palermo (1975); S. Mandelstam, Phys. Rept. 23, 245 (1976). 

[10] Z. F. Ezawa and A. Iwazaki, Phys. Rev. D 25 (1982) 2681; T. Suzuki and 
I. Yotsuyanagi, Phys. Rev. D42, 4257 (1990). 

[11] For a review see: M. N. Chernodub and M. I. Polikarpov, arXiv:hep-th/9710205. 

[12] V. G. Bornyakov et al, Nucl. Phys. Proc. Suppl. 106 (2002) 634. 

[13] V. G. Bornyakov et al, arXiv:hep-lat/0310011. 

[14] J. D. Stack, S. D. Neiman and R. J. Wensley, Phys. Rev. D 50 (1994) 3399. 
[15] H. Shiba and T. Suzuki, Phys. Lett. B 333, 461 (1994). 

[16] G. S. Bali, V. G. Bornyakov, M. Muller-Preussker and K. Schilling, Phys. Rev. D 54 
(1996) 2863. 

[17] S. Ejiri, S. Kitahara, T. Suzuki and K. Yasuta, Phys. Lett. B 400 (1997) 163. 

[18] Y. Mori et al, Nucl. Phys. A 721 (2003) 930; V.G. Bornyakov et al, 
arXiv:hep-lat/0301002; Nucl. Phys. Proc. Suppl. 119 (2003) 703. 

[19] B. Sheikholeslami and R. Wohlert, Nucl. Phys. B 259 (1985) 572. 

[20] S. Booth et al. [QCDSF-UKQCD collaboration], Phys. Lett. B 519 (2001) 229; M. 
Gockeler et al., to be published. 

[21] C. R. Allton et al. [UKQCD Collaboration], Phys. Rev. D 65 (2002) 054502. 

[22] R. G. Edwards and U. M. Heller, Phys. Lett. B 462 (1999) 132. 

[23] K. Jansen and R. Sommer [ALPHA collaboration], Nucl. Phys. B 530 (1998) 185. 

[24] A. Hasenfratz and F. Knechtli, Phys. Rev. D 64 (2001) 034504. 

[25] V. G. Bornyakov, D. A. Komarov and M. I. Polikarpov, Phys. Lett. B 497 (2001) 
151. 

[26] F. Brandstater, U.-J. Wiese and G. Schierholz, Phys. Lett. B 272 (1991) 319. 



-27- 



[27] T. A. DeGrand and D. Toussaint, Phys. Rev. D 22, 2478 (1980). 

[28] J. Smit and A. van der Sijs, Nucl. Phys. B355 (1991) 603. 

[29] T. Suzuki, S. Ilyar, Y. Matsubara, T. Okude and K. Yotsuji, Phys. Lett. B347 
(1995) 375; Erratum: ibid. B351 (1995) 603. 

[30] H. Satz, Nucl. Phys. A 681 (2001) 3. 

[31] M. N. Chernodub, E.-M. Ilgenfritz and A. Schiller, Phys. Lett. B 547, 269 (2002). 
[32] S. Aoki et al. [CP-PACS Collaboration], Nucl. Phys. Proc. Suppl. 73 (1999) 216. 
[33] H. Satz, arXiv:hep-ph/0111265. 

[34] F. Gliozzi and P. Provero, Nucl. Phys. B 556, 76 (1999). 

[35] T. Suzuki and M. N. Chernodub, arXiv: hep-lat/0207018, hep-lat/0211026. 

[36] S. Digal, P. Petreczky and H. Satz, Phys. Lett. B 514 (2001) 57. 

[37] S. Kratochvila and P. de Forcrand, Nucl. Phys. B 671 (2003) 103. 

[38] C. DeTar, O. Kaczmarek, F. Karsch and E. Laermann, Phys. Rev. D 59 (1999) 
031501. 

[39] L. D. McLerran and B. Svetitsky, Phys. Rev. D 24 (1981) 450. 
[40] S. Nadkarni, Phys. Rev. D 33 (1986) 3738. 

[41] O. Kaczmarek, S. Ejiri, F. Karsch, E. Laermann and F. Zantow, hep-lat/0312015. 
[42] M. Liischer, P. Weisz, JHEP 0207 (2002) 049. 

[43] M. Gao, Phys. Rev. D 40 (1989) 2708; see also: Ph. de Forcrand, G. Schierholz, H. 
Schneider and M. Teper, Phys. Lett. B 160 (1985) 137. 

[44] G. S. Bali et al. [TXL Collaboration], Phys. Rev. D 62 (2000) 054503. 

[45] C. Y. Wong, Phys. Rev. C 65 (2002) 034902. 

[46] A. Hasenfratz, R. Hoffmann and F. Knechtli, Nucl. Phys. Proc. Suppl. 106 (2002) 
418. 

[47] O. Kaczmarek, F. Karsch, E. Laermann and M. Lutgemeier, Phys. Rev. D 62 (2000) 
034021. 

[48] G.S. Bali, K. Schilling and C. Schlichter, Phys. Rev. D51 (1995) 5165; 

R.W. Haymaker, V. Singh, Y.-C. Peng and J. Wosiek, Phys. Rev. D53 (1996) 389; 

[49] V. Singh, D. A. Browne and R. W. Haymaker, Phys. Lett. B 306 (1993) 115; 

M. Zach, M. Faber, W. Kainz and P. Skala, Phys. Lett. B358 (1995) 325; G.S. Bali, 
Ch. Schlichter and K. Schilling Prog. Theor. Phys. Suppl. 131 (1998) 645. 



- 28 - 



[50] A. Hasenfratz, Phys. Lett. B 476 (2000) 188. 



[51] A. Hart and M. Teper, Phys. Lett. B 371 (1996) 261; 

V. G. Bornyakov and G. Schierholz, Phys. Lett. B 384 (1996) 190; 
M. N. Chernodub and F. V. Gubarev, JETP Lett. 62 (1995) 100. 

[52] S. Hioki, S. Kitahara, S. Kiura, Y. Matsubara, O. Miyamura, S. Ohno and 
T. Suzuki, Phys. Lett. B 272 (1991) 326 [Erratum-ibid. B 281 (1992) 416]. 

[53] Y. Nakamura et al, arXiv:hep-lat/0309144. 



- 29 - 



